home *** CD-ROM | disk | FTP | other *** search
/ SGI Hot Mix 17 / Hot Mix 17.iso / HM17_SGI / research / lib / trace.pro < prev    next >
Text File  |  1997-07-08  |  2KB  |  72 lines
  1. ; $Id: trace.pro,v 1.5 1997/01/15 03:11:50 ali Exp $
  2. ;
  3. ; Copyright (c) 1996-1997, Research Systems, Inc.  All rights reserved.
  4. ;       Unauthorized reproduction prohibited.
  5. ;+
  6. ; NAME:
  7. ;       TRACE
  8. ;
  9. ; PURPOSE:
  10. ;       This function computes the trace of an N by N array.
  11. ;
  12. ; CATEGORY:
  13. ;       Linear Algebra.
  14. ;
  15. ; CALLING SEQUENCE:
  16. ;       Result = TRACE(A)
  17. ;
  18. ; INPUTS:
  19. ;       A:      An N by N real or complex array.
  20. ;
  21. ; KEYWORD PARAMETERS:
  22. ;       DOUBLE: If set to a non-zero value, computations are done in
  23. ;               double precision arithmetic.
  24. ;
  25. ; EXAMPLE:
  26. ;       Define an array, A.
  27. ;         A = [[ 2.0,  1.0,  1.0, 1.5], $
  28. ;              [ 4.0, -6.0,  0.0, 0.0], $
  29. ;              [-2.0,  7.0,  2.0, 2.5], $
  30. ;              [ 1.0,  0.5,  0.0, 5.0]]
  31. ; 
  32. ;       Compute the trace of A.
  33. ;         Result = TRACE(A)
  34. ;
  35. ;       The result should be: 3.00000
  36. ;
  37. ; REFERENCE:
  38. ;       ADVANCED ENGINEERING MATHEMATICS (seventh edition)
  39. ;       Erwin Kreyszig
  40. ;       ISBN 0-471-55380-8
  41. ;
  42. ; MODIFICATION HISTORY:
  43. ;       Written by:  GGS, RSI, January 1996
  44. ;-
  45.  
  46. FUNCTION Trace, X, Double = Double
  47.  
  48.   ON_ERROR, 2
  49.  
  50.   Sx = SIZE(X) 
  51.   if Sx[0] ne 2 then $
  52.     MESSAGE, "Input array must be 2-dimensional."
  53.  
  54.   if Sx[1] ne Sx[2] then $
  55.     MESSAGE, "Input array must be square."
  56.  
  57.   if N_ELEMENTS(Double) eq 0 then $
  58.     Double = (Sx[Sx[0]+1] eq 5) or $
  59.              (Sx[Sx[0]+1] eq 9)
  60.  
  61.   ;TOTAL(DoubleData, Double = 0) returns a double-precision result. Cast the
  62.   ;result to COMPLEX or FLOAT if Double = 0.
  63.   if Double eq 0 and Sx[Sx[0]+1] eq 9 then RETURN, $
  64.     COMPLEX(TOTAL(X[LINDGEN(Sx[1]) * (Sx[1]+1)], Double = Double)) else $
  65.   if Double eq 0 and Sx[Sx[0]+1] eq 5 then RETURN, $
  66.     FLOAT(TOTAL(X[LINDGEN(Sx[1]) * (Sx[1]+1)], Double = Double)) $
  67.   else RETURN, $
  68.     TOTAL(X[LINDGEN(Sx[1]) * (Sx[1]+1)], Double = Double)
  69.  
  70. end
  71.  
  72.